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Abstract- This paper presents a finite element model using cube elements to characterize electromagnetic field. The 
modeling approach uses an optimized solution via edge element method. This work shows that the computational 
cost is significantly reduced. The computational aspect is guaranteed by GNU Octave, a high-level interpreted 
language primarily intended for numerical computations. GNU Octave provides capabilities for the numerical 
solution of linear and nonlinear problems. It also provides extensive graphics capabilities for data visualization and 
manipulation. 
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I. INTRODUCTION 

A variety of numerical methods have been proposed for computing the electromagnetic properties in frequency 
and time domain as finite element method (FEM), finite difference (FD), method moment (MM), boundary integral 
method (BIM). The FEM remains the most used one and this is due to its ability to incorporate the different types of 
boundary conditions and to model the inhomogeneous material with complicated geometry and properties [1, 2, 3, 

4]. 

This work uses GNU Octave [5] software to perform an efficient electromagnetic study using cubic finite 
element [6,7]. The use of edge elements is an efficient method way to solve this kind of problem [8, 9, 10]. Octave 
is a high-level interpreted language, essentially, intended for numerical computations. It provides capabilities for the 
numerical solution of linear and nonlinear problems. Octave provides a rich graphics capabilities for data 
visualization. It is used through an interactive command line interface (CFI), and it can also be used as programming 
language using the non-interactive programs. The Octave aspect is quite similar to Matlab so that most programs are 
easily portable. 

Octave is freely redistributable software, used for intensive numerical computation with an interactive 
environment that is compatible with Matlab. It is customizable with user-defined functions written in Octave’s 
language or C++, C, Fortran or others [13]. 

A waveguide [11, 12] is an electromagnetic transmission line used in microwave communications, 
broadcasting, and radar installations. A waveguide consists of a rectangular or cylindrical metal tube or pipe. The 
electromagnetic field propagates lengthwise. 

An electromagnetic field (E,H) can propagate along a waveguide in various ways. Two common modes are 
known as transverse-magnetic (TM) and transverse-electric (TE). This paper demonstrates the efficiency of cubic 
element to analyze such structures using minimum of elements and computational resources via GNU Octave free. 

The remainder of this paper is organized as follows. In Section 2, we introduce the 3D Maxwell problem mixed 
FEM. In Section 3, we give the cubic formulation of the posed problem. In section 4 some numerical results and 
finally the main conclusions in last section. 
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II. 3D MAXWELL PROBLEM 




Figure 1: Transmission line containing three domains. The length, 1=1 5cm 

We consider the 3D problem of a waveguide transmission line. It consists of two parallel conducting plates of 
width a separated by height b by a dielectric material containing three domains {_£ n ,p n ),n = 1,2,3. Example of such 
line a can be used as a microstrip line in the context of microwave integrated circuits. 

For a full wave analysis we deal with the Maxwell equations in the frequency domain: 



V ' x E = sB 


(1) 


VxH = J + sD 


(2) 



Where s = ja), E the electric field intensity, H the magnetic field intensity, B the magnetic field flux, D the electric 
field flux and J, the is the current. 

The constitutive relations are given by: 





B = pH 


(3) 




D = eE 


(4) 


The Ohm law gives: 








J = gE 


(5) 



From (1-5), we establish the vector wave equation in terms of the electric field E: 

V X p- 1 X E + s(<r + 5 s)E = ~sj imp (6) 

A computation using the finite element method is performed in a finite. In order to truncate the volume of the 
computational domain the Silver-Miiller condition is applied as an absorbing boundary condition. It is given by: 

n x V x E = jkE tan 

where k is the wave number in free space and E tan is the tangential electric field on the outer boundary surface. 

This ABC preserves the sparsity and symmetric features of the final matrix. It is exact for normal incidence. 

Let denote Q the computational domain and T the outer boundary. As usual with FEM, we define the space of the 
work: 

v = {ue (L 2 (n)) 3 , v xu e (L 2 (n)) 3 } 

A weak formulation of the problem is obtained after multiplying the vector wave equation by a test function F in V: 
(s(a + se)E, F) n +(^ _1 V x E, V x F) n +(s A /£M =T £’- F) r = 0 (9) 
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Where On denotes the scalar product in V. 

To solve the system (9) numerically, the domain is discretized with cubic elements. The electric field can be written 
in terms of basis functions associated with the edges of these elements: 

E = Ze u Nlj 

Using test function F the same as interpolation functions (Galerkin method) and supposing that a = 0, we get the 
following system: 



(Gj£ — 0 ; 

< (Gj = M + 5 § + 5^ M, 

, e = ( e ij) ’ 

Where e^ j , is defined as circulation of the electric field along the edge which is the path between two nodes 

i and j: 

e i,i = f { i, n Edl ’ 

The matrices are defined by: 

— ^i,j> Njc,v}si 

%JV fe ,n)J = r 

%(ij),(fc,n)} = (M -1 V X Nfj, V X Nl n ) a 
In the above, matrices entry comes from edge (i,_/) and edge (k, h). 

HI. CUBIC FEM FORMULATION 

Let us Consider a cubic element given in Fig. 2, whose side length is L and whose center is at (xc, yc). 




Figure 2: Oriented cubic element, 6 faces 12 edges and 8 nodes. 



The cubic edge element can be written as follow [1]: 



Nh = 



L 3 



\y c +■ 



-’][ 



z r +- - Z\ 



X, 
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w « = Z3l 



-z r + - + z 



NU = 



Nh=^Ur + 



L 3 



y c + 2 -y]b . 

♦he*§- 



X, 



y< 



v 2,3 



L 3 
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IV. NUMERICAL EXAMPLE 

The transmission line the length 1=1 5cm is excited by a plane wave E(z)=\.e~ kz y on Z = 0 surface, where k is 
the wave number, and z the direction of wave propagation. The frequency used in this study is y=2.45 GHz. The 
finite element mesh for this example was generated by our own 3D mesh generator. We note that 



If we take the permittivity and the permeability of the studied domain e 1 = £ 2 = £ 3 = 1 and p 1 = p 2 = l 1 3 = 1 (the 
vacuum), and l x = l 2 — l 3 = 5 cm. The figure 3 represents the mesh of the studied which contains 244 nodes, 484 
edges and 60 elements. Fig. 4 and Fig. 5 represent modulus of electric field and the magnetic field respectively for 
as function of z. The obtained results are compared with analytical ones. As we can see from these figures, we obtain 
the same results for different value of z. 
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Figure 3: The mesh for waveguide, with 244 nodes, 484 edges and 60 elements. 



The Magnetic Field 




z[m] 



Modulus of E[V/m] 




Figure 4: The modulus of electric field. 



The Magnetic Field 





Figure 5: The modulus of magnetic field. 
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V. SPARSITY OF MATRICES 



A. Elementary matrices 

The elementary matrices are 12 by 12 matrices, because of 12 nodes in the cube. 
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Figure 6: Sparsity pattern of the elementary matrices 



B . Global matrices 

The global matrices are really sparse. This, sparsity corresponds to systems which are loosely coupled. The Sparse 
matrices have significant advantages in terms of computational efficiency. Using sparse matrices to store data can 
both save a significant amount of memory and speed up the processing of that data. 
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Figure 7: Sparsity pattern of the global matrices: M, § and M 
VI. CONCLUSION 

The cubic elements discretize the Maxwell wave equations. The solution is efficient an accurate. The problem is 
optimized using spare matrices. GNU Octave gives a good tool to study a complexe FEM problem, and offer good 
alternative of open source software. 
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APPENDIX 

A. The mesh- How-to 

Technically, the cube element is represented in (Fig. 10). 



A(ij+1,k+l) A(i+1 j+1 ,k+1 ) 




The matrix A, allows managing the nodes of the cube element. The following script shows how to store node 
number and the position of each node in the cube, and finally in the whole mesh, where D x and D y are the 
discretisation 

■ Nodes 

iNode=l; 

for ix=l:Nx+l 

for iy=l:Ny+l 
A(ix,iy,l)= iNode; 

Node.x(iNode) = (ix-l)*Dx; 

Node.y (iNode) = (iy-l)*Dy; 

Node.z(iNode) = Z(l); 
iNode = iNode +1; 
end; 

end; 

Node.Total = length(Node.x); 

■ Edges 

%othey-edges 

iEdge=l; 
for ix=l:Nx+l 

for iy=l:Ny 
k=k+l; 

Edge.NodeNum ( iEdge, 1) = A(ix,iy,l ); 

Edge.NodeNum ( iEdge, 2)= A(ix,iy,l)+1 ; 
iEdge = iEdge+1; 
end; 

end; 

%the x-edges 
for iy=l:Ny+l 

for ix=l:Nx 

Edge.NodeNum ( iEdge, 1)= A(ix,iy,l ); 

Edge.NodeNum(iEdge,2)= A(ix+l,iy,l); 
iEdge = iEdge+1; 
end; 

end; 

Edge.Total = iEdge- 1; 
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■ Elements 

The manage a connectivity between element and nodes, we use the matrix j ^ 

PNT=[]; 
for kx=l :Nx 
for ky=l:Ny 
K=[kx,ky]; 

PNT=[PNT;[K, k] ]; 
end; 
end; 



for iElt=l:Nx*Ny 
i=PNT(iElt,l); 
j=PNT ( Elt,2); 
k=PNT ( Elt, 3 ); 
Elt.NodeNum(iElt,l) 
Elt.NodeNum(iElt,2) 
Elt.NodeNum(iElt,3) 
Elt.NodeNum(iElt,4) 
Elt.NodeNum(iElt,5) 
Elt.NodeNum(iElt, 6) 
Elt.NodeNum(iElt, 7) 
Elt.NodeNum(iElt,8) 
Elt.NumDom(iElt) 
end; 



= A(i,j,k); 

= A(i+l,j,k ); 

= A(i+lj+l,k); 

= A(ij+l,k); 

= A(i,j,k+1); 

= A(i+l,j,k+l); 

= A(i+lj+l,k+l); 
= A(i,j+l,k+l); 

= 1 ; 



B. Analytic solution 

The analytical solution E=(Ex,Ey,Ez) of such configuration is given by: 

r, r \ sin (JkL-kz) 

E y (z) = sinOcL) ’ 

From the equation (1) we can compute H={Hx,Hy,Hz) using Maxwell Faraday equation. 

^ (Z) = sgsin ( fc L) C ° S(feL “ fcz) / sin m, 
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